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C/2 , We study the shear-induced dilTusion effect and the transition to irreversibility in suspensions 

, under oscillatory shear flow by performing an analysis of the entropy production associated to the 

^ , motion of the particles. We show that the Onsager coupling between different contributions to the 

, entropy production is responsible for the scaling of the mean square displacement on particle diam- 

C , eter and applied strain. We also show that the shear-induced effective diffusion coefficient depends 

Q_) , on the volume fraction and use Lattice-Boltzmann simulations to characterize the effect through 

the power spectrum of particle positions for different Reynolds numbers and volume fractions. Our 
study gives a thermodynamic explanation of the the transition to irreversibility through a pertinent 



(N 

> , 

^\ analysis of the second law of thermodynamics 

(N ■ 
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00 I. INTRODUCTION 

o ■ 

. ^ When a suspension of non-Brownian particles is subjected to an oscillatory shear flow, the dynamics of the 
' particles presents a transition to irreversibility which has been recently observed in experiments [l|, 3, Q . In 

5t 1 these experiments, the suspension of polymethylmethacrylate (PMMA) particles having sufflciently large sizes 
(diameter d ~ 230^m) , is contained in a cylindrical Couette cell and taken out of equilibrium by applying an 
oscillating shear flow proportional to 7 cos{ujt), where 7 — uj^q with 70 the applied strain and lo the characteristic 
frequency of the oscillation. At small enough Reynolds numbers it is observed that the motion of the particles 
is oscillatory and reversible, according to a classical result of hydrodynamics j^. When increasing the Reynolds 
number or the concentration of particles, the trajectories of the particles become chaotic and then their reversible 
behavior is lost. This effect is manifested through a shear-induced diffusion which has been characterized through 
the mean square displacement (MSD) of the particles 

ABB. 

The MSD scales in the form: (Ax^) ^ (f^t, 
thus allows one to deflne an effective diffusivity scaling like D ^ (P'j. 

Characterizing the motion of the particles through the MSD clearly suggests that a statistical description of 
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their dynamics is possible. Previously, this description was offered in Ref. [5] by postulating a diffusion equation 
in which the diffusivities have been constructed by analyzing the temporal behavior of the position correlation 
function of the particles. This approach allows the use of direct experimental measurements or simulation 
results in order to describe particular systems ^j. Other theoretical and numerical studies have characterized 
the relation between the transition to a chaotic motion of the particles with the shear-induced diffusion effect 

In this article, we offer a general description of this shear-induced diffusion effect and the associated transition 
to irreversibility which is based on the application of the second law of thermodynamics, and on previous works 
devoted to analyze the dynamics of a suspension of Brownian particles in the presence of flows [7, 8]. 
We calculate the entropyproduction of the system in the phase space of the particles and find the corresponding 
Onsager couplings 0,3,01. One of these coupHngs is responsible for the dependence of the diffusion tensor on the 
imposed velocity gradient, even in the limit of small specific thermal energy: ksT /m 0, with /c^T the thermal 
energy and m the mass of a particle Q . This dependence of the diffusion tensor on the velocity gradient leads 
to the shear-induced diffusion effect that depends crucially on hydrodynamic interactions lOl. Illl. Il 
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and the breaking of the fluctuation-dissipation relation (FDR) at mesoscopic level, 13 liaM 1111181,11 
The article is organized as follows. In Sec.|ll]we use non-equilibrium thermodynamics in phase space to formulate 
the mesoscopic model based on a Fokker-Planck equation. Sec. Ill is devoted to deriving a Smoluchowski 
equation having an effective diffusion tensor accounting for the shear-induced diffusion observed in experiments. 
In section IV we present Lattice-Boltzmann simulations characterizing the shear-induced diffusion effect via 
the power spectrum of particle movements affected by hydrodynamic interactions. Finally, in Sec. V we discuss 
our main results. 



II. MESOSCOPIC ENTROPY PRODUCTION FOR THE DYNAMICS OF A SUSPENSION IN 

EXTERNAL FLOW 

We consider a suspension of N non-interacting spherical particles of radius a and mass m in a fluid which moves 
with velocity u"(r, t). Since the system is in contact with a heat bath that evolves in time, it is necessary 
to determine the physical nature of the coupling forces in order to adequately describe its dynamics. This 
objective can be achieved by taking into account two factors. The first one is that the evolution of the system 
can be described at mesoscopic level by means of the normalized A^-particle probability distribution function 
P'^^^(r^, i), that depends on the instantaneous positions {r}^ = (ri, r/v) of the particles and their velocities 
{u}^ = (ui, ...,un) through the phase space vector = {{f}^ , {u}^). The second factor takes into account 
the fact that the interactions between the system and the heat bath involve dissipation. This suggests the use 
of the nonequilibrium entropy s{t) as a thermodynamic potential from which the entropy production a(t) can 
be calculated, and used to obtain the explicit expressions for the coupling forces 23[. 
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To proceed in systematical way, we will first notice that the probability distribution function satisfies the 
conservation law 



i=l i—1 

where \7p. represents the gradient operator with respect to the position vector and J^. is a diffusion current 
defined in phase space. Integration of Eq. ([T]) over the phase space coordinates , under the assumption that 
P^^) and Js- vanish at the boundaries, leads to the continuity equation: dp/dt —V - {pv) in which the average 
density field p{r,t) of the suspended particles is defined by 14 1 

N 

p(f, t)^m J2 P^''Hr'',mn - r)dr^, (2) 

and the mean velocity field v{r, t) is 

N 

v{r,t) = -m / y w,pW(r^,i)^(r. - r)dT^ . (3) 
Here dV^ ~ d{u\^ d{r}^ is the volume element in the phase space of the particles. 

One of the purposes of this section is to derive explicit expressions for the currents Ja^ which at this point 
implicitly contain the mentioned coupling forces between system and bath. Once these expressions are obtained, 
the evolution equation for p(^) can be written. As we have mentioned, Jp can be obtained from the entropy 
production of the system which follows from the Gibbs entropy postulate |l4j | 



N 



J2pm\n^6{n-f)dr^, (4) 

where 5s is the entropy change with respect to a local equilibrium reference state characterized by the local 
equilibrium distribution function 



l.eq. ^ 



(5) 



Here /^s is the local equilibrium chemical potential per mass unit and v^ ^ = v^{ri,t). 

Following the rules of mesoscopic nonequilibrium thermodynamics [ol, we take the time derivative of Eq. ^ 
and use H]), an integration by parts over F^-space assuming that the fluxes vanish at the boundaries, leads to a 
balance equation for the entropy s in which the entropy production contains three contributions: a = X]j=i ^j - 
The first contribution is related to the diffusion process in {ti}-space 

-1 = -?/E4-^'5(.-.-r)dF-, (6) 

i — 1 



4 



where the noncquihbrium chemical potential /i(r^, t) is given by 

N 



kRT 



(7) 



The second contribution comes from diffusion of particles with respect to the mean velocity, with diffusion 
current Ji — {ui — Vi)P^^'': 

N 

(8) 



/iV 

i—1 



The third contribution corresponds to diffusion with respect to the flow velocity whose current is Jp = (ui 

N 



(9) 



where Fi = dvf /dt is a non-stationary force related with the variation of the fluid velocity with time. 
According to the second law of thermodynamics, the entropy production of the system must be positive definite 
(7 > for irreversible process. To satis fy th is condition, nonequilibrium thermodynamics establishes linear 
relationships between currents and forces 23|. In particular, for Jg. we obtain 



N 



N 



N 



where the tensors dij, eij and Cj, are related to the Onsager coefficients L 



and L,, in the form 



(10) 

Q 

(11) 



The Onsager coefficients obey Onsager's relations in which time-reversal symmetry must also be applied to the 
external drive: Lu^rj — —Lnuj, 2M- From Eq. pUj) it follows that the system of particles is coupled to the heat 
bath by means of thermal and entropic forces (first term on the right hand side of the equation) and mechanical 
forces (last two terms of the equation). 

Substituting Eq. (|10p into the continuity equation for the probability ([1]), we arrive at the multivariate Fokker- 
Planck equation describing the evolution of the iV-particle distribution function 

QpiN) N N 



.EV.,-(u.P(-))=E4 



dt 



m ouj J 



-1 



(12) 



where we have used Eqs. ([7]) and ((TO)) assuming that the coefficients are symmetric tensors. Finally we introduced 
the friction tensor /3j ■ leading to the relation Q 



(13) 
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It is important to mention that the combination eij ■ "S/p-ifj, entering in Eqs. ()12|) and (I13p . imphes that the 
fluctuation-dissipation theorem (FDT) connecting the drift and diffusion terms of the Fokker-Planck equation 
is no longer valid due to the presence of the shear 0, lisl . [l6, 22|. This important consequence following from 



Eq. is related to the shear-induced diffusion effect, as we will show in the next section. In the case of a 
diluted suspension, similar results for the diffusion term of the generalized Fokker-Planck equation have been 
obtained by means of the kinetic theory of gases in Ref. [21 1. 

The coefficient eij is related to the force exerted on the surface of a particle moving through a fluid under flow 
conditions. For a spherical particle, e^j has been calculated explicitly in terms of the generalized Faxen theorem 
in Ref. |25|, and used in Ref. [8,] to obtain tij = eoc^j with 

eo = ^^a'/3o'(l + 2aa), (14) 

where /3o — drnja/m is the Stokes friction coefficient per mass unit and t] the viscosity of the fluid, a = 
{—iio jv)^!"^ is the inverse viscous penetration length, u: is the frequency and v the corresponding kinematic 
viscosity 13] ■ The tensor e,j is related to the friction tensor /Jj^ and obeys the relation ta = I with I the 
unit tensor. In Ref. 8] it has been shown that is related to inertial effects due to the change in time of vq 
and has the form — C^Sij, with C — Pp/pf, Pp the density of the particle and pf the density of the heat 
bath. From Eq. (fT4| it follows that ep incorporates finite-size effects on the dynamics of the system through the 
surface term and frequency-dependent corrections to the diffusion coefficient through a. 
Eqs. (fT2|l - (fT4|) imply that the diffusion coefficient in velocity space, [ksT /ra) dij, does not vanishes in the limit 
ksT/m and, therefore, also imply that the non-thermal contribution to the diffusion coefficient associated 
with the Onsager coefficient e.y may have important consequences on the dynamics of a non-Brownian suspension 
of particles. We will show these consequences in the next section. 

The friction tensors /3y are affected by the hydrodynamic interactions among particles and their dependence 
can be inferred from its relation with the mobility tensors /i^ : /S^ • p.^j = iSij. At lower-order approximation, 
the multipole expansion of /Ij ■ takes the form 



III 



(15) 



Here fij and f,j-^ are the unit relative vectors between particles and between particle j and the wall, rij = — f*,- 
is distance between particles whereas r^^ is the magnitude of the vector that points from sphere i to the mirror 
image with respect to a wall of sphere j. For i — j, Eq. (fT5|) reduces to well-known results for the mobility of a 
particle in the presence of a wall: p = (3^^ (I — Bia/l), with I its distance to the wall. The coefficient Bi may 
take different values depending on the direction of the motion of the particle with respect to the plane of the 
waU [26 1 . 
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III. SHEAR-INDUCED DIFFUSION 



We will analyze in this section the diffusion regime occurring at times t ^ Pq^- Since experiments and 
simulations give the self-diffusion coefficient [if [23I , we will focus our description on the dynamics of a single 
particle whose reduced distribution function pk{f,t) — m /P*^^-'(5(ffc — rjdT^ satisfies the continuity equation 



dpk 
dt 



-V • (pkVk), 



(16) 



which follows by integrating Eq. (|12|1 over the phase space of the particles and where Vk (r, t) — 
mp^^ JukP^^'>S{rk~r}dT^. 

The Smoluchowski equation for pk can be derived after calculating the explicit expression for pkVk- This task 
can be carried out by calculating the evolution equations for the momentum field pkVk and for the pressure 
tensor of the fc-th particle, defined as 

Pfe(r,t) = m J{uk - Vk){uk - Vk)P'^''^5{rk - r)dr^. (17) 

Following the method indicated in Ref . Q , we take the time derivative of the definition of pk {r, t) , use Eq. 

in the result and perform an integration by parts assuming that the currents vanish at the boundaries. After 

rearranging terms we arrive at the equation 



Pk 



(18) 



where we have used the expression C^^j = C,15ij and defined the convective derivative dk/dt = d/dt+Vk-d/df, and 
the two-particle reduced distribution function p^'^\rk,ri,t) = m J P^^^dT^^'"^. Here, c?r^^^ the phase-space 
volume element of the N — 2 particles including the measure associated to the velocities Uk and Ui. In 
the superscript indicates a dependence on fk and r*; . The right hand side of Eq. (jlSp represents the total drag 
force exerted by the fluid on the particles. The first term on the right hand side is the friction force including 
the presence of hydrodynamic interactions that modifies the local value of the velocity field in terms of the 
distribution of particles in the system. 

The evolution equation for the pressure tensor can be derived by following a similar procedure: 



dt 



/3ol + Vi7fc + iv-Vfel 



2kBT 



PoPk 



a^(3o{ l + 2aa) ViT^ 



(19) 



N 



i—l,i^k 



where to obtain the second term on the right hand side, we have used the expressions (jl3p and (|14p for i ~ j. The 
upper symbol s means the symmetric part of a tensor. The last term on the right hand side of this equation 
contains the contribution of the hydrodynamic interactions to the pressure tensor of the particle k. This 
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contribution enters through the cross-correlation functions of the local velocities indicating how hydrodynamic 
interactions modify the stresses in the system. In order to obtain a closed expression for the pressure tensor 
Pfc , it thus becomes necessary to calculate the evolution equation for the cross-correlation function 



N 



(20) 



The evolution equation for C^i can be obtained in a similar way as we did to derive Eq. (fT9| . The result is 

N 



dt 



(21) 



N 



E 4- • (^i ' - Wi - ^)P^'^'5(r, - r)drudn 

where k ^ i. To derive Eqs. (|19p and (|2ip . we have neglected the contributions arising from higher order 
moments of the time-ordered hierarchy since they relax faster than the ones present in those equations Q . 
At times t ^ we can obtain from (flQl) the following constitutive equation for Pfc 



m 



-Pk 



1 



a^/3o(l + 2aa) V^jjj 



(22) 



AT 



/^<7'y E /3..-(4'^-^)(^^f -t^')P^'^'5(rl-^1drl.drl 



where we have used that {Vvk)ijf3Q ^ 1, V • Vkf^Q ^ 1, in accordance with the experiments 
For times i > /3o'\ Eq. ^ can be rewritten in similar form as Eq. (|22p by extracting the term j = k from 
the sum on the second term at the right-hand side, multiplying the resulting relation by (3q^ (3^^ and performing 
the sum over i. One obtains the expression 



m 



N 

E Pki ■ dkip''^^S{fk - r)drkdfi (23) 

N 

E • (4'^ - ^f'k){€^ - ^)p^'^S{rk - r^dvkdn, 

i—l,i^k 

where we have used Pf.k — /^olj kept terms of the order (a/rij)'^ and neglected terms of the order (a/rij)^ and 
higher . This approximation is valid up to intermediate volume fractions of the suspended particles. 



A. Effective medium approximation 

The Smoluchowski equation for pk can be obtained from (fT6|) . (fT8|) . ((22|) and ([23| by assuming an effective 
medium approximation in which the test particle k performs its motion in a fluid incorporating the effects 



of hydrodynamic interactions in average form . In our description, this assumption does not consider the 
possibiUty of direct collisions among particles, and thus we expect that it is valid up to intermediate volume 
fractions of the suspended particles. Operationally, this approximation can be implemented by substituting 
the averages appearing on the right-hand side of equations HH]), and ([^5]) by integrals over a continuum 
variable. Thus, the last term on the right-hand side of Eqs. l(22l) and ((23)) may be written as 

E ■ (4'^ - «1)(4'^ - ^)P^'^m r)drkdn J /3(r') • C (f - r' , t)dr^ , (24) 

i—l^i^k 

where we have introduced the non-local velocity cross-correlation function C (f — r' ,t). In similar form, for the 
momentum field we have the relation 

E ■ (^'^ - ^.)P^^^S{rk - r^drkdn ^ j ■ {v - ^).., p^^) (F_ t)df , (25) 

i—l^i^k 

where we have defined (w — i}*) -.^, = v{r — r' ,t) ~ 'if'{r — r' ,t). Using Eqs. and (HH), the right-hand side of 
Eq. ([23]) can be rewritten as 

E ■ ^Lp^^^Sivk - r)drkdn ^ j d*{^)p^^\r~ r', i)dr' = (26) 

l(r') •|(f')p(2)(f_^,t)dr'- -^a'/?2(lH-2aa)[Vi/']t • / [/f(r') • l(r')]tp(2) _ ^-'^ i)^^^ 

where the upper f means the transpose of a tensor and we have used the fact that ViP does not depends on 
position. The tensor a ir') — (3{r') ■ air') has been defined to simplify the notation in subsequent relations. 
In the long-time limit the substitution of Eqs. ([24| and (l26l) into (122]) yields the constitutive equation for the 
pressure tensor 



(27) 



where we have assumed that the second term on the left hand side of Eq. ((2T|l may be neglected, see Eq. ([22|l . 
Using now Eq. (j25p in (|18p and taking the long-time limit t S> Pq^' from the momentum equation we obtain 

V • I ~ - / l(r') • {v- i/') , t)dr^ + p^F. (28) 



The explicit expression for the constitutive relation of the diffusion current pv of the particles follows by assuming 
that the spatial variation of the velocity field is small enough in order to make the expansion 

v{r- r', t) - v°{r~ r*, t) ~ v{f) - v°{f) + O [(Vi^)^j , (29) 

where we have used ~ 1?° + 0(V In p), Q. 
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Factorizing the two-particle distribution function in the form 29]: p^'^'^ (f— r' , t) ~ p{'r)g(f—r',t), with g^r — r' ,t) 
the two-particle correlation function, we define the effective quantities 



Bir, t- 4>) = j /3(r') • /3(r ')g(r - r', t- 4>, T)dr' and E{r, t; 0) = /Jq-^ J Pir') ■ i{r')g{f- , t; cb, T)dr', 

(30) 

where we have taken into account the fact that the two-particle correlation function g may in general depend 
on the volume fraction and the temperature, 3ol |. 

After using these results into Eqs. (|^5|) - (P5|) . we obtain the following constitutive equation for the diffusion 
current 



pv ~ pv 



m 



\7-A]p + pCB • F - D{r, t) ■ Vp, 



where we have identified the effective diffusion tensor 

D(r, t) ^ kBTimp H (1 + 2a a) 

6 



Here p = Pq fl is the effective mobility tensor and we have introduced the dimensionless tensors 



1 and A 



m 



-y(3o{l + 2aa) iE-l)-Vv 



(31) 



(32) 



(33) 



Eqs. ([50)1 - (|32p show that the transport coefficients in the coarse-grained description contain hydrodynamic 
interactions in effective form through the configurationally averaged tensors B and E. When Vw* — 0, then 
the diffusion tensor reduces to: D{r, t) = {ksT /•m)p{r, t; 0), that is, in equilibrium the diffusion tensor depends 
on hydrodynamic interactions through the mobility tensor p, as expected. The coarse-graining performed in 
this section thus leads to incorporate the contribution of hydrodynamic interactions in the dissipation of the 
reduced system by modifying the diffusion coefficient making it anisotropic and position and time dependent. 
Substituting now Eq. (j3ip into (jl6p . we finally obtain the Smoluchowski equation 

dp 



dt 



pvo - p^a' 1 + M • ./ - PCf^o' 1 + A • ^ 



V • [D-Wp 



(34) 



where we have defined the force / due to hydrodynamic interactions as / = — (fcsr/m)V • A. 
Eqs. (l32l) - ((34)) constitute the main result of this section. The effective diffusion tensor D contains two contribu- 
tions. The first one depends on the thermal energy per mass unit (ksT/m) and is therefore related to Brownian 
motion whereas the second one does not depends on thermal fluctuations. 
In the limit of vanishing specific thermal energy, ksT/m 0, Eq. ([32]) leads to 

„2 



D 



6 



(1 + 2aa) 



1 + p] ■ {E -I) -V^ 



(35) 



This expression has the same scaling on particle diameter and shear rate as that observed in experiments and 



Stokesian dynamics simulations [ll. l27j. The presence of /i(r, t; cj)) and E{r, t; (j)) implies that the shear-induced 
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diffusion is mediated by hydrodynamic interactions. We then conclude that hydrodynamic interactions are 
responsible for the randomization in the motion of the suspended particles when an oscillatory strain is applied 
on the system. 



B. The mean square displacement 

We will assume that the fluid velocity vo(f,t) = f- j{t) is imposed along the x direction with j{t) the time- 
dependent shear rate whose only non- vanishing component is jyx = 'ycoslujt). The shear rate is related to 
the applied strain 70 by 7 = 70^^^. For convenience, we will assume that inertial effects are negligible and that 
effective mobility fl and E do not depend on time and position. This hypothesis is valid when the distribution 
of the suspended particles does not changes significatively, that is when g{f' — r, t; (p) ^ g(r '; </>). 
The MSD of particle's position vector can be calculated by taking the time derivative of the expression 

{r') = / {x'+y')pdr, (36) 



where we have considered the two-dimensional case. Substitution of Eq. (|34p into the result and an integration 
by parts leads to 

^(r2) = 2jcosiL0t) {xy)it) + 2Tr[B], (37) 

where {xy){t) = Jxypdf. In similar form, we must derive the evolution equations for {xy){t), {x'^){t) and 
{y'^){t). After solving the obtained set of differential equations, for low shear rates and frequencies (7 < 1, 
Lu < 1) we may expand the MSD in a power series of 7 and lu to obtain 



(r^) ~ ADf) + jlyy] t + ^jlxy[E - l) (f jt 



24 



(38) 



where Dq = ksT/mfSo is the one-particle diffusion coefficient and d ~ 2a having assumed p,xy — Pyx and 
E = El. In the limit ksT /m — s- when the particles are non-Brownian we obtain 

{r^)^\p.ymE{^)-Ad^t. (39) 

When Eq. ([5^ is expressed in terms of the number of oscillations n of the imposed flow, with t = 27rn/cj, 
it gives: (r^) ~ ^PxyiE — l)d'^jon. This relation shows that hydrodynamic interaction introduce a volume 
fraction dependence of the MSD of the particle, thus giving an explanation based on thermodynamic arguments 
for the scaling relation obtained in the experiments [l| . 

In the limit of massive particles, from Eq. ((34| it is also possible to derive the evolution equation for the average 
position of the particle defined through R{t) = J rpdf. By taking the time derivative of this definition and 
integrating by parts one obtains 

j^R{t) = m ■ Vvo{t) - CPo^6,{R-4>) ■ j^^Mt) + yf^2(^;<^) • V«o(0, (40) 
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where we have introduced the quantities Gi = ((^1 + ' ^^^d G2 = (V • [(^1 + ■ (-E" — 1)]) and the 
bracket indicates an average over p. Eq. (|40p is a nonUnear equation for R{t) in which the nonhnearities are 
a consequence of hydrodynamic interactions through the terms Gi and G2. In a first approximation, the last 
term at the right hand side of Eq. (|40p estabhshes that hydrodynamic interactions become significant when the 
Reynolds number defined by Re = dv^jv satisfies the condition 

Re > h{(j>), (41) 

where v is the kinematic viscosity of the heat bath and the function h{(f)) takes into account that (G2)y is a 
function of cj). In obtaining this relation we have neglected the second term at the right hand side of ([3D]) since 

is a very small quantity and scaled time with uj and lengths with d. 
Eqs. ([55)1 and ([1T|) show that the transition to the irreversibility is mediated by hydrodynamic interactions that 
introduce a dependence on volume fraction of the shear-induced diffusion coefficient. 

IV. SHEAR-INDUCED DIFFUSION FROM LATTICE-BOLTZMANN SIMULATIONS 

In this section, we analyze the shear-induced diffusion effect by means of Lattice-Boltzmann simulations. This 
method allows us to study the dependence of the effective diffusion coefficient as a function of the relevant 
parameters of the problem, the Reynolds number Re and the volume fraction (f). The power spectrum of the 
components of the trajectories of the particles is used to show that the randomization of particle movements is 
due to an increasing number of modes produced by hydrodynamic interactions as Re and 4> increase. 
We use the two-dimensional model D2Q9 for the Lattice-Boltzmann method with the BGK approximation 
[3^. In this model, the space is discretized in a two dimensional square lattice with nine velocities (c^, « = . 
allowed. The particle distribution functions /(r, i), at site r and time t evolve according to the equation 



31 



1 



/.(r + cO-/.(r,i)^-- f,ir,t)~ft'\r,t) , (42) 

T L J 

where t is the dimensionless relaxation time related to viscosity and /j^'^'^'' are the local equilibrium distribution 
functions, 

l-HSc, •u+^(c, •u)2-^w2 (43) 



fi = WiP 



In this equation, Wi — 4/9, 1/9, 1/36 are the weights associated to the lattice 33| for each set of velocities 
|ci| = 0, 1, \/2 and p and u are the density and velocity defined by 

p(r,<)=^/,(r,i), u(r,t) = i^/,(r,t)c,. (44) 

i ^ i 

The viscosity is related to the dimensionless relaxation time by = Cg(T — 1/2), where = l/\/3 is the speed 
of sound in the D2Q9 model. 
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iqo V--^- I ^ ^ ^ ^ ^ ^ ^ 1 
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Re 

FIG. 1: Effective diffusion coefficient as a function of Re (symbols) for a fixed particle concentration of = 0.14 and 
/* = 10.0. The solid line is a linear fit with slope of 5.68 x 10~^. 



The no-slip boundary conditions are simulated on the solid particles and the torques and forces are also evaluated 
to update the particles position at all times [3J|. The interactions among particles are implemented with the 



method proposed in Ref. [35] and with the corrections proposed in Ref. 34|. The walls of the cavity use the 
bounce-back boundary condition, which consists in reversing the incoming particle distribution function after 
the stream process. 

The numerical simulations are carried out in a cavity of H* = 11.33 and W* = 44.66 where the dimensions 
were scaled with the radius of the particle. The relaxation time and the radius of the particles are kept constant 
in all simulation at t — 20.0, r — 4.5 and /* — 10.0 and we have varied the Reynolds number and the volume 
fraction. The dimensionless frequency was scaled with the magnitude of the shear rate. 

In the first set of numerical simulations, we fixed the particle concentration (f> = 0.14 and the dimensionless 
frequency /* = 10, varied the Reynolds number and determined the effective diffusion coefficient D as shown 
in Fig. [TJ A linear dependence of D^ff on Re has been obtained to a good approximation, as expected from 
experiments and theoretical results. Using the definition Re — (Pj/v, from Eq. (|39p we may also obtain the 
linear relation for the effective diffusivity as a function of the Reynolds number: D^f / = flxy ((p) — \]vRe/(!>. 
The data obtained from simulations can be used to give a rough estimate of the magnitude of hydrodynamic 
interactions in terms of the parameter ilxy{4>)[E{(f)) — 1]. Given d — 9 and ly = 11.25 we obtain p,xyi(t')[E'{(t>) — 
1] ^ 10^"*. This small value indicates that 3D Lattice-Boltzmann simulations are required in order to do a 
quantitative comparison with experiments. 

In order to better discern the mechanisms leading to the shear-induced diffusion, we calculated the power 
spectrum (PS) of the components of the trajectories for fixed (j) and Re, as it can be respectively seen from 
Figs. [2] and [3 

In Fig. [2] (a), we present the PS of x{t) for four different values of Re and (f> = 0.14. The insets represent 
three different trajectories for the same particle with the same initial condition for Re = 0.01, Re = 0.07 and 
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Re — 0.08. For lowest Reynolds number (solid line), the PS presents more pronounced peaks located at the 
excitation frequency and its harmonics. The trajectory in the inset (i) (solid line) corresponds to this spectrum 
and shows a very regular behavior. From the PS for the case of Re = 0.04 (dotted line) it follows that the 
dynamics in the x-component keeps the main peak and harmonics at the same position than in the previous case, 
but small peaks start to appear between the harmonics. This is a consequence of the hydrodynamic interactions 
between the particles that give rise to new frequencies in the dynamics of the system. In the case Re — 0.07 
(dashed line) there is a shift of the harmonics and the new frequencies are better defined. The corresponding 
trajectory (dashed line) is shown in the inset (ii). Finally, for Re — 0.08 (short-dashed line) the harmonics 
disappear and the energy is more homogeneously distributed for frequencies larger than the one imposed. The 
corresponding trajectory (dotted line), shown in the inset (Hi), is irregular. In Fig. [2] (b) we present the PS 
corresponding to the y— movement of the particle for different Re at a fixed (j). From this set of PS we can 
appreciate that both the x and movements are coupled with the exciting frequency, as well as the fact that 
the harmonics have a small shift and the energy is distributed in more frequencies. The PS of x{t) and y{t) 
shown in Figs. [21(a) and (b) agree with theoretical results in two ways. First, they indicate a coupling between 
different modes, like in Eq. (|37p . Second, the power spectra reflects that hydrodynamic interactions become 
important only for Reynolds numbers larger than a certain values, as established by Eq. (|4T|) . 
The theoretical results ([55)1 and ([55)1 indicate that the shear-induced effect is caused by hydrodynamic in- 
teractions and, although it is not explicitly shown, one then expects that the effect also depends on particle 
concentration. This was confirmed by performing another set of simulations keeping Re = 0.08 constant and 
varying (j) to analyze the influence of hydrodynamic interactions on the dynamics of the particles through the 
PS of x{t) and y{t). The results are shown in Figs. [3lja) and (b), respectively. For the lower volume fraction 
(f) = 0.00636 (one particle) the PS for x{t), Fig. ^a), has only one peak at the exciting frequency. This result 
is expected from theory, Eq. (|40p . since for sufficiently small strains and frequencies the leading contribution 
corresponds to the first term on the right hand side, which is proportional to the applied strain. For the y(t), 
Fig. [Sljb) , the peak is shifted to the right from the exciting frequency, showing a weak coupling between the x and 
y motions. For (j) ~ 0.01270 (two particles), x{t) presents harmonics which disappear at frequencies much larger 
than the exciting one whereas for y(t) the peak at the exciting frequency appears with an incipient presence of 
harmonics, implying that hydrodynamic interactions are weak. For (j) — 0.02540 (four particles) three harmonics 
can be identified for both x{t) and y{t). For cj) — 0.05090 (eight particles ) a larger number of harmonics can 
be identified in x{t) with about half of the energy in comparison with the exciting frequency. For the y(t), 
the harmonics are clearly identified and have the same energy as the exciting frequency. New frequencies arise 
between the harmonics implying that hydrodynamic interactions introduce new modes in the dynamics of the 
particles. This result is expected from the evolution equation (j40p which predicts that for strains or particle 
concentrations larger than a certain critical value, new modes will appear in the behavior of R(t) due to the 
contributions of the nonlinear terms. Finally, for <j) = 0.102 (16 particles), the only peak for the x{t) and y{t) is 
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FIG. 2: Power spectrum for the same particle at different Reynolds numbers for (a) the x-movement and (b) the y- 
movement. The cross in the insets indicates the initial position of the particle, the same in the three cases. Notice that 
maximum displacement increases with increasing Reynolds number. 



located at the exciting frequency and the energy is now homogeneously distributed in more frequencies. These 
results indicate that an increase in particle concentration enhances the effects of hydrodynamic interactions 
which in turn are responsible for distributing the energy in a growing number of modes. They also reflect the 
fact that hydrodynamic interactions also introduce a dependence on volume fraction through the dependence 
on the volume fraction of the coeflacients jlxy{4>) and E{(j)). 
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Power spectrum for difFerent values of (j) at Re — 0.08 for a) the x and b) j/— movement of the trajectory of a particle. 



V. CONCLUSIONS 



In this paper, we have analyzed the shear-induced diffusion effect in suspensions under oscillatory shear by 
means of a thermokinetic theory based on the calculation of the entropy production at mesoscopic level. We have 
found that an Onsager coupling between thermal and non-thermal effects containing hydrodynamic interactions 
is responsible for this effect. 

By calculating the entropy production of the iV-particle system in contact with a heat bath and identifying the 
corresponding forces and currents, we have derived a Fokker-Planck equation for the iV-particle phase space 
distribution function having an effective diffusion coefficient in which the coupling between thermal and non- 
thermal effects breaks the validity of the fluctuation-dissipation theorem. After contracting the description over 
velocity space, assuming that the diffusion regime is well established and using an effective medium approxi- 
mation, we derived a Smoluchowski equation for the single particle distribution function containing an effective 
diffusion tensor incorporating the cross-effects of the imposed flow and hydrodynamic interactions, Eq. (|32|) . In 
the limit ksT/m — > 0, this diffusion tensor yields the same scaling on particle diameter and applied strain as 
formd in experiments and simulations, D ~ d? flxy{4>)[E{(l)) — l]jo/6. Our analysis may also explain the influence 
thermal noise on this effect. 
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The Lattice-Boltzmann simulations we performed were used to show in more detail the dependence of the 
shear-induced diffusion on the Reynolds number and the particle volume fraction. To this effect we have used 
the power spectrum of the components of particle trajectories. As expected from theory, we found a linear 
dependence of the diffusion coefRcient Z) as a function of the Reynolds number Re. When increasing Re for 
fixed (f) as well as cj) for fixed Re, we obtained that the power spectrum shows an increasing contribution of new 
modes, period doubling and finally, for large Re and (j), a loss of characteristic frequencies indicating a stochastic 
behavior induced by hydrodynamic interactions, in accordance with theoretical results expressed through the 
mean square displacement (|39|) . It is important to emphasize that the 2D simulations performed clearly indicate 
that up to intermediate volume fractions the hydrodynamic interactions are responsible for the shear-induced 
diffusion effect, in agreement with the theoretical prediction. A more precise description could be carried out 
from 3D simulations in order to perform a quantitative comparison with experiments. 

The expression for the shear-induced diffusion coefficient in Eq. (I39p and the condition ([H]) show that the 
transition to the irreversibility is due to hydrodynamic interactions and that it depends on the volume fraction, 
in accordance with experiments and theoretical simulations. Therefore, our study gives a theoretical explanation 
of the shear induced effect and the transition to the irreversibility based on the pertinent analysis of the entropy 
production at mesoscopic level. 
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